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ABSTRACT 

We generate high quality grids interactively on a computer workstation through opti- 
mization for some simple yet commonly encountered geometric shapes. 


1. INTRODUCTION 

Algebraic and elliptic grid generation are the two most widely used methods for gen- 
erating grids. Algebraic grid generation is known for fast interactive grid generation. The 
grid smoothness and lack of skewness, however, are usually not assured. On the other 
hand, elliptic grid generation can generate smooth grids but the time it takes to generate 
them is usually too long for interactive response. It is naturally desirable to combine the 
advantages of both methods and to generate high quality grids (ref. 1) fast and efficiently 
on a computer workstation. This requires people to strike a proper balance between the 
grid generation speed and the grid quality. We have pursued this idea for some simple 
yet commonly encountered geometric shapes. Many of them have a symmetry that can be 
used. We have been able to obtain high quality grids fast and interactively. Our method 
is to analyze each case carefully and to take full advantage of the region symmetry (if any) 
so that optimization is needed only for a part of geometric region under consideration. 
We have combined various grid generation methods, such as a fast elliptic grid genera- 
tion technique (combining the control point form of algebraic grid generation with elliptic 
grid generation, ref. 2) and a transfinite interpolation with orthogonal boundary terms, to 
achieve an optimization for each individual case. After such an optimization, we have been 
able to generate excellent single- and multi- block grids quickly. This achievement enables 
us to generate certain single- and multi-block grids interactively with reasonable response 
even on modestly powered workstations. 
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2. GRID WITHIN A RHIOMBUS 


A simple transfinite interpolation (TFI) will not generate an orthogonal grid in general. 
A more general TFI with boundary derivative terms can improve grid orthogonality. In 
certain cases, grid orthogonality can be achieved along edges of a face grid or on surfaces 
of a volume grid. This, in general, requires a careful analysis for each geometry region 
encountered. In this section, we present an example of carrying out such a detailed analysis. 


We consider a surface grid bounded within a quadrilateral. Since each of four boundary 
edges is a straight line, a bicubically blended interpolation formula (i.e., with derivative 
information) can be used, 


r (u,v) = B T (u) S B(v) , 

where S is a 4 x 4 matrix 


(i) 


/ r(0,0) r(0, 1) 

I r(l,0) r(l, 1) 

r«(0,0) r u (0, 1) 

\r u (l,0) r„(l,l) 


r„(0,0) r v (0,l) \ 

r„(l,0) r B (l,l) 
r«w(0,0) r uv (0, 1) I 
r„„(l,0) r u „(l, 1) / 


and B_ is a 4-component vector 


/#i(«)\ 

B 2 (u) 

Bz(u) 

\B a {u)/ 


/I — 3u 2 -(- 2u 3 \ 

| 3 u 2 — 2 u 3 I 

I u — 2 u 2 + u 3 I 
V — u 2 + u 3 ) 


(2) 


( 3 ) 


In Eq. (1), u and v are parametric variables. The parametric region of 0 < u < 1, 
0 < v < 1 corresponds to the area within the quadrilateral. In Eq. (2), r(0,0), r(l,0), 
r(0, 1), and r(l,l) are four corners of the quadrilateral (i.e., the surface grid), four r„ s 
are partial derivatives with respect to the variable u at four corners, four r„’s are partial 
derivatives with respect to the variable v at four corners, and r uv ’s are second-order mixed 
derivatives (also called twists) at four corners. In Eq. (3), the blending functions Bj{u) 
have the properties: B\{u) + B 2 {u) = 1, Bj( 0) = 8j \ , Bj( 1) = 8j 2 , dBj(0)/du = Sj 2 , 
dBj(l)/du = Sj 4 - 


In grid generation, the positions of four corners r(0,0), r(l,0), r(0,l), and r(l,l) 
are known. However, the values of 8 first-order derivatives and 4 second-order derivatives 
are not given. With a proper choice of these derivatives, one can generate a better grid 
than those generated by the simple TFI or, for a quadrilateral, by the simple bilinear 
interpolation. The idea is to require that the grid direction at the boundary edges point 
to a desired direction. For example, on the v = 0 edge, the grid direction going inside the 
quadrilateral can be easily found from Eqs. (l)-(3) to be 


r v (u,0) = Bi(u)r v (0, 0) + B 2 (u)r v (l, 0) + R 3 («)r„„(0,0) + B 4 {u) r ut ,(l,0) 


( 4 ) 
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Uy requiring it to be perpendicular to the v = 0 edge (a straight line here), it is possible 
to achieve boundary orthogonality at certain places on the v = 0 boundary edge. 

In the following, we take a rhombus as an example to illustrate how a detailed analysis 
can dramatically improve the grid orthogonality. A rhombus is a special parallelogram with 
opposite sides and opposite angles equal but having no right angle. For a rhombus centered 
.at origin (on the z = 0 plane), its four corner positions can be written as (see fig. 1) 

r(l,l) = — r(0,0) = h j , r(l,0) = -r(0,l) = w i . (5) 



Figure 1. — A rhombus centered at 

origin. 

From symmetry considerations, we know (i) that the magnitudes of r u (0, 0), r„(0, 0), 
t*( 1,1), and r w (l,l) must be equal, and (ii) that the magnitudes of r«(l,0), r„(l,0), 
a*„(0, 1), and r„(0,l) must be equal too. In other words, we have the following relations: 

r u (0,0) = r u (l, 1) = au , 

(6a) 

r„(0,0) = r„(l, 1) = av , 

(6b) 

r u (l,0) = r u (0, 1) =fiu , 

(6c) 

r v (l, 0) = r„(0, 1) = /3v , 

(6d) 


where 
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u = (un + hj)/y/w 2 + h 2 , v — (~wi + hj)/\/w 2 + h 2 


( 7 ) 


are the unit vectors along the directions of r(l,0) — r(0,0) and of r(0,l) — r(0,0), respec- 
tively. It is clear that both a and (3 are positive: a > 0 and (3 > 0. Also, because of 
symmetry, (i) twists r(0,0) and r(l,l) do not have x component whereas twists r(l,0) 
and r(0, 1) do not have y component, and (ii) the magnitudes of the twists at each pair 
of opposite corners are equal but point to opposite directions. Namely, twists satisfy the 
relations 


r uv (l,l) = — r uv (0, 0) = a j , 


(8a) 


r* v (l,0) = — r uw (0,l) = b i . (8b) 

Expressions (6) and (8) are the most general form allowed (by the symmetry requirements) 
for a rhombus. In other words, the symmetry consideration reduces the number of inde- 
pendent first-order derivatives from 8 to 2, and also reduces the number of independent 
second-order derivatives from 4 to 2. Now we only need to decide 4 scalar quantities: a, 
ft, a, and b. For the rhombus shown in fig. 1, Eq. (4) reduces to 


r v (u,0) = [aBi(u) + /3i? 2 (u)]v - aBz(u) j + bB^(u) i . (9) 


Without loss of generality, we assume in the following discussion that h > w. When 
h > w, the angle between the two unit vectors u and v, which is the internal angle at 
corner r(0,0), is less than 90°, since v • u = ( h 2 — w 2 )/(h 2 + u; 2 ) > 0. To determine 
the 4 quantities a, /?, a, and 6, we first consider an (unnormalized) orthogonality function 
f(u) = r v (u, 0)-u, which is a function (third-order polynomial) of the parametric variable u . 
Its value at the two ends r(0, 0) and r(l, 0) of the v = 0 edge is positive: /( 0) = av • u > 0 
and /( 1) = (3ir • u > 0. We want to force the grid to be orthogonal at the middle of 
the v = 0 edge (but not overly doing this enforcement). Specifically, we require that the 
orthogonality function f(u) reaches its minimum value at u = | and that the minimum 
value of the orthogonality function is zero (for the interval 0 < u < 1). These two 
requirements lead to two conditions for the orthogonality function /(w), 

df( i) 

*r =0 ’ «*>-"• < l °) 

From these two conditions, we can express the twists a and 6 in terms of the first-order 
derivatives a and /?, 


(5a — /3)(h 2 — w 2 ) 
h\/w 2 + h? 

(lla) 

(5 P - a)(h 2 - w 2 ) 
w\/w 2 + h 2 

(lib) 
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Eliminating the twists a and b in the orthogonality function /(«) with the help of Eqs. 
(H)» /(«) becomes 

/(») = iv • u [(a + /3)(2u - 1)’ - (a - «(2« - l) 5 ] , (12) 

which indeed reaches its minimum value of 0 at v = | (for the interval 0 < u < 1 ) for any 
positive a and /?. 

Next, we ask r„ and r u to be orthogonal at the center of the rhombus r(|, |) = 0. 
Due to symmetry, the angle between r„(|, |) and j is always equal to the angle between 
r»(|, |) and j. Thus, asking r„ and r„ to be orthogonal at the center of the rhombus is 
equivalent to requiring that r„(|,|)‘j = r„(|,|)* (— i) and 

r«(|,|)-i = r u (i,i).j. (13) 

Equation (13) leads to a third relation among the derivatives a, /?, a, and b, 


r»(M)-(j-i) = 


h — w 


4 y/w 2 + h? 


(6\/ w 2 + h 2 — a — /3) + 


a + b ^ 3((3 — a)(w + h ) 


16 


8\/iu 2 + h 2 


(14) 


Substituting Eqs. (11) into Eq. (14), we can eliminate a and b and obtain an expression 
relating a to j3, 


(w + 5 h)/3 + 24 wh(h — w)(u/ 2 + h 2 ) l ! 2 

a = 

5w + h 

Now we can still choose the value of fi as it fits. There is, however, a restriction. In the 
case of h > w under discussion, the twists a and b , as determined by Eqs. (11), should not 
be negative. Otherwise, the resulting grid will not look nice. 



In practice, we find that a vanishing twist at the pair of obtuse angles of the rhombus, 


i.e. 


b = 0 , 


is a good choice. This leads to 


(16) 


a = 5/3 

according to Eq. (lib). It follows from Eqs. (15) and (17) that 

„ h(h — w) 

P — . = , h > w . 

Vw 2 + h 2 

Substituting Eqs. (17) and (18) into Eq. (11a), we find that correspondingly 

24 (h — w)(h 2 — w 2 ) 


a = 


w 


2 + h 2 


(17) 


(18) 


(19) 
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Equations (16)-(19) thus uniquely determine all 12 derivatives in Eq. (1). 

Figure 2 shows three rhombuses inside a regular hexagon of side length 2w. Two 
rhombuses with a pair of edges parallel to the x-axis can be obtained from the third 
rhombus (oriented as in fig. 1) with a simple rotation of il20 around the corner r(l,0) 
on z = 0 plane. For each rhombus in fig. 2, h = y/3w and Eqs. (18) and (19) give 
/? = 0.634iu and a = 8.785w. 



Figure 2. — A three-block grid inside a regular hexagon. 


3. OTHER EXAMPLES 

When the geometric region for a surface grid is no longer a quadrilateral (as in most 
cases), a TFI formula with boundary information terms, i.e., a bicubically blended TFI 
(see, for example, ref. 3), should be used to replace Eq. (1). By a careful analysis for 
each geometric region or shape under study, one can properly choose the directions and 
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magnitudes of boundary derivative terms and thus dramatically improve the grid quality 
and smoothness. In this way, a surface grid with much improved grid orthogonality can 
be generated. For a volume grid, one can also include boundary derivative terms into the 
often-used volume grid FTI formula and thus generate a volume grid with improved grid 
orthogonality. We have carried out this approach for some simple yet commonly encoun- 
tered geometric shapes besides the rhombus case discussed above in section 2. Besides this 
pure algebraic grid generation approach, we have also combined the control point form of 
algebraic grid generation with elliptic grid generation (ref. 2). Here we first elliptically 
generate a sparse net of control points and then use the control point form of algebraic grid 
generation to generate a much denser grid. As a result, we generate a smooth grid with 
less computational time and storage space. This approach also enables us to obtain a high 
quality grid in an interactive environment on a computer workstation. We are incorporat- 
ing these new grid generation techniques into a software product GridPro™/IGgen from 
Program Development Corporation, so that a user will be able to readily generate high 
quality grids interactively on a computer workstation. Currently, a collection of boundary 
geometries are provided and a user can choose among them. In addition, the user can also 
specify the grid type wanted, such as O-type, C-type, H-type, and/or multi-block grids. 
The user can also specify the grid density and clustering pattern. Once the user makes 
a choice, our software immediately generates high quality grids. In this section we show 
some of the grids generated interactively using GridPro/IGgen. 

Figure 3 shows a single-block 2D grid for a passage which has a circular extra space. 
Due to symmetry, we only need to carry out grid orthogonality optimization for a half, a 
quarter, or even one eighth of the grid. 



Figure 3.— A single-block 2D grid inside a 2D pipe which has symmetric circular extra 
space. 
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Figure 4 displays an H-type 2D grid passing around a cylinder whose cross section 
looks like a super circle. There are two blocks. One block is the mirror image of the other 
block about the horizontal symmetry line. For each block, there exists a vertical symmetry 
line. So we only need to do optimization for one quarter of the whole grid. 



Figure 4.— A 2-block grid passing around a supercircle-like object. 
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Figure 5 shows a single-block volume grid inside one half of a torus. Five surfaces 
(except the top surface) are displayed. Symmetry of the geometry shape enables us to 
focus on one sheet of the volume grid and do grid orthogonality well. Other sheets can be 
obtained by rotations of different angles. 



Figure 5.- — Half of a single-block volume grid inside a torus. 
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Figure 6 shows two pictures of a 7-block volume grid inside a sphere. Because of 
symmetry, we only need to perform optimization for one of six side blcks and the center 

block. 



Figure 6. — Two views of a 7-block volume grid inside a sphere. 
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These examples clearly show the smoothness and lack of grid skewness that can be 
reliably expected in our new interactive environment. 


4. SUMMARY 

We have combined several new grid generation methods to obtain high quality grids 
in an interactive environment. Grids can be generated reliably, and grid quality is high in 
general. Currently, we are developing a grid generation program called GridPro T M /IGgen 
to take advantage of these new advances in grid generation. Our goal is to generate a range 
of high-quality single- and/or multi-block grids interactively on a computer workstation. 
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